Two-phase stretching of molecular chains 
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While stretching of most polymer chains leads to rather featureless force-extension diagrams, 
some, notably DNA, exhibit non-trivial behavior with a distinct plateau region. Here we propose a 
unified theory that connects force-extension characteristics of the polymer chain with the convexity 
properties of the energy vs. stretch profile of its individual monomer subunits. Namely, if the 
inter-atomic interaction potentials within the monomer are such that its is described by a uniformly 
convex (concave down) energy function, the stretching of the entire chain is single phase and uniform. 
However, if the monomer stretching energy as a function of extension has a non-convex (concave up) 
region, the stretching of the polymer chain separates into two phases: weakly and strongly stretched 
monomers. We use simplified planar and 3D polymer models to illustrate the basic principles of the 
proposed model. Specifically, we show that when the secondary structure of a polymer is mostly due 
to weak non-covalent interactions, the stretching is two-phase, and the force-stretching diagram has 
the characteristic plateau. This scenario is realized in a-helix, consistent with earlier observations 
based on all-atom molecular dynamics simulations. Using a coarse-grained model of the DNA with 
realistic potentials, we show explicitly that the double helix exhibits the two-phase behaviour under 
stretching, in agreement with single molecule experiments. The predicted end points of the force- 
extension plateau are in close quantitative agreemnt with experiment, while a reasonable qualitative 
agreement is obtained for the tension value at the plateau. Detailed analysis demonstrtates that 
breaking of Watson-Crick bonds is not necessary for the existence of the plateau, although some of 
the bonds do break as the double-helix extends. 



I. INTRODUCTION 



Early on, the mechanical properties of DNA attracted 
the interest of both physicists and biologists, because of 
their importance to numerous biological processes, such 
as DNA transcription, gene expression and regulation, 
and DNA replication. Besides, mechanical properties 
are essential while constructing nanodevices containing 
DNA. 

In the past decade new methods appeared which made 
it possible to stretch single molecule and measure a 
stretching force with high precision. In these experiments 
one end of the molecule is attached to the surface, an- 
other one is attached to a latex bead which is pulled by 
optical or magnetic tweezers [1-4]. The main result of 
such experiments is a force-extension diagram. Complex 
behavior has been revealed by elasticity studies of indi- 
vidual double-stranded DNA molecules. Under normal 
conditions several different regimes have been observed 
[5]. Until distance between molecule ends is less than 0.8 
of its contour length the force (chain tension) increases 
linearly with extension increasing, according to Hooke's 
law. It should be mentioned that rigidity in this regime 
is rather small, 2-^3 pN force is sufficient to stretch a 
molecule for 90% of its contour length [1]. The force in 
this regime has purely entropic nature. When the end- 
to-end distance is close to the contour length the force 
increases significantly [1, 2]. This effect is explained by 
intrinsic elasticity, i.e. related with structure deforma- 
tion, which begins to contribute significantly when the 
contour length of the molecule is approached. Here elas- 
tic response is due to both entropic factor and intrinsic 



elasticity. So when end-to-end distance is less than con- 
tour length molecule behaviour is described by ordinary 
persistent model [5]. When the molecule is subjected to 
forces of 65-^120 pN, it suddenly changes form, stretching 
up to 70% beyond its canonical B-form contour length. 
Various models for the structure of this so-called S-form 
DNA have been proposed, but its nature remains ob- 
scured. When higher forces are applied molecule divides 
into two strands and afterwards covalent bonds break, 
chemical structure is disrupted. 

Single molecule stretching experiments have been car- 
ried out also for polypeptide molecules. [6-8]. The com- 
mon polypeptide has complex secondary structure, and 
its force-extension diagram is rather complex; but similar 
diagrams with typical force-extension plateau have been 
observed experimentally for synthetical alpha-helices [7] 
and myosin molecules [6]. 

Force-extension plateau observed in single DNA 
molecule experiments is usually explained by either grad- 
ual unzipping of the double helix (force-induced ten- 
sion) [3, 4, 9, 10], or cooperative transition of the whole 
structure into a new form called S-form where hydro- 
gen bonds between complementary bases remain unbro- 
ken. The latter is supported by the fact that transition 
occurs in a narrow range of forces (about several pN) 
and is reversible [1, 2], as well as by existing of the sec- 
ond transition under significantly higher forces (about 
150 -T- 300 pN), which is attributed to non-equilibrium 
melting [11]. The version of the force- induced tension is 
supported firstly by estimates of free energy which reveal 
that under room temperature the melted state is more 
favorable due to entropic effects [9], secondly by coin- 
cidence of force-extension diagrams for double-stranded 
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and single-stranded DNA on certain range [9] , thirdly by 
experiments which visualize double- and single-stranded 
segments using specific fluorescent labels [12]. Computer 
simulations also don't enable to determine the nature of 
the force-extension plateau: different investigation sup- 
port existing of S-form DNA, where base pairing remains 
intact [1, 13, 14], as well as possibility of force-induced 
melting [15, 16]. 

In the case of polypeptide force-extension plateau is 
attributed to the alpha-helix unwinding [6, 17, 18]. each 
site is supposed to have two stable states which have 
different free energy — winded and unwinded. Then this 
representation can be rather easily described by Izing 
type models [19, 20]. This phenomenon is called helix- 
coil transition. 

Assumption that two phase (strongly and slightly 
stretched) occur in polymer molecule stretching process 
has been suggested previously for DNA molecule [1, 19- 
21]. Usually the strongly stretched state was supposed to 
consist of melted base pairs, where hydrogen bonds are 
ruptures [9, 19, 20]. But a model has been supposed [21] 
where hydrogen bonds remain in strongly stretched state, 
which is stable and corresponds to an energy minimum. 

The paper considers mechanical stretching of three dif- 
ferent molecular chains: plane zigzag chain, helical chain 
and DNA double helix in coarse-grained representation. 
Our results reveal that for linear polymers non-uniform 
stretching can take place only when non-valent interac- 
tions contribute significantly to its secondary structure 
stabilization. It is the case for polypeptide a-helix which 
is stabilized by hydrogen bonds and for DNA double he- 
lix which is stabilized by stacking interactions and hydro- 
gen bonds between complementary nitrogen bases. Thus 
plateau in a force-extension diagram is a consequence of 
relation between non-valent and valent interactions. 

It should be noted that current paper consider only 
mechanical energy (instead of free energy which is con- 
sidered in most papers). However even consideration 
of this constituent only enables to explain a typical 
force-extension diagram. Besides, computer simulation 
of DNA double helix dynamics under the temperature 
300 K, performed within the framework of coarse-grained 
representation, indicates that main result remains correct 
independently of thermal fluctuations. In the stretched 
state hydrogen bonds are weakened significantly, so cer- 
tain hydrogen bonds may rupture due to thermal fluctu- 
ations. However these hydrogen bonds ruptures can not 
be considered as a cause for the plateau on the force- 
extension diagram. 

Assumption about two-phase stretching have been pro- 
posed previously [6, 17, 21], but the stretched state was 
always considered as stable, that is corresponding to the 
free energy minimum. The current paper demonstrates 
that such an assumption is not necessary for two-phase 
scenario of stretching (the right end of the plateau may 
be not the minimum of the energy function) . 




Figure 1: Dependence of uniformly stretched molecular chain 
energy E (curve 1) from relative longitudinal step extension 
Ap. The straight line (line 2), connecting points Api and Ap2 
and being the tangent in these points, describes a convex hull 
of the function E(Ap), Apo stands for the inflection point. 



II. THE GENERAL MECHANISM 



Aforementioned features of a force-extension diagram 
for DNA molecule and polypeptide alpha-helix can be 
explained in a generalized way. The chain structure is 
stabilized by interactions of two types: valent and non- 
valent. The former describe bond, angle and torsion 
deformations. The latter describe interactions between 
atoms which are not connected by rigid valent interac- 
tions. These interactions include electrostatic and van 
der Waals interactions, hydrogen bonds in polypeptide 
alpha-helix and stacking interactions between neighbor- 
ing base pairs in DNA. The main feature of the non- 
valent interactions potentials V(r) is existing of the in- 
flection point which is a left end of an interval where 
deformation energy increasing slows significantly. If non- 
valent interactions contribute significantly to uniformly 
stretched chain energy E(Ap) (Ap > is relative ex- 
tension), then the function E(Ap) also may have an 
inflection point Apo. In this case the function E(Ap) 
is no longer convex downwards; its convex hull can be 
obtained by connecting two points (Api, E(Api)) and 
(Ap 2 , E(Ap 2 ))) (0 < Api < Ap < Ap 2 ) by a common 
tangent - see Fig. 1. This line describes more favorable 
way of stretching when one part of the chain is in strongly 
stretched state with extension Api, and another one is 
in slightly extended state with the extension Ap 2 - While 
this non-uniform stretching chain extension is due to in- 
creasing of the strongly stretched sites, so the tension 
remains constant. Such non- uniform stretching scenario 
explains completely the existing of the force-extension 
plateau. 
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The energy of non-valent interactions between second 
neighbors is given by the potential 



n-2 n—1 n n+1 n+2 n+3 

Figure 2: Scheme of the zigzag molecular chain. 



III. STRETCHING OF THE 2-D ZIGZAG 
MOLECULAR CHAIN 

Let's explain the mechanism of the non- homogeneous 
two-phase stretching on example of 2D zigzag chain de- 
picted on the Fig. 2. Let a molecular chain have a form 
of plane zigzag with fixed x-step h x and y-step h y . All 
atoms of the chain are situated in the xy-plane. In the 
equilibrium state n-th atom has coordinates x n = nh x , 
y n = [1 + (—l) n }h y /2. A zigzag is described by the dis- 
tance between neighboring sites pa = (h 2 + hy) 1 / 2 and 
zigzag angle 9$ [x-step is h x ~ posin(9o/2), y-step is 
h y = p a cos(0o/2)]. 

The 2D zigzag energy is given by the following hamil- 
tonian: 

H = + + + U ( n) + W(r n )}, (1) 



with M being an atom mass, (x n ,y n ) - n-th atom coor- 
dinates. The potential V(p n ) describes the energy of in- 
teraction between neighboring n-th and (n + l)-th atoms 
depending from the distance p n = (a^ + b^) 1 ' 2 , where 
a„ = x n+ i - x n , b n = y n+1 - y n . This potential corre- 
sponds to the deformation energy of a valent bond be- 
tween neighboring atoms. The potential U(9 n ) describes 
the deformation energy of the angle between three atoms 
(n — l)-th, n-th and (n + l)-th, which is given by 

n = arccos[-(a„a„_i + 6 ra 6„-i)/p„-ip„]. 

The last term W(r n ) describes weak non-valent interac- 
tion between n-th and (n + 2)-th atoms (second neigh- 
bors), with r„ = [(x n+2 - x n ) 2 + (y n+2 - y n ) 2 } 1/2 being 
a distance between these atoms. 

Let's consider dimensionless zigzag model where equi- 
librium bond length is po = 1, zigzag angle is 9q = 
arccos(-l/3) = 109.47° (so that x-step is h x = 0.8165, 
y-step is h y = 0.5774). 

Rigid bond deformation is given by 



V(p) = ±K(p-p Q ) 2 , 

with bond rigidity being K = 2. 
Angle deformation is given by 



(2) 



U{9) = ei (cos 9 -cos 6» ) 2 , (3) 
where t\ > is angle deformation coefficient. 



W(r) = £ 



rp - d 
r — d 



- 1 



(4) 



where eo > is interaction energy, tq = 2h x = 1.633 
is equilibrium length, d = 0.5 is a diameter of the inner 
hard core. For the sake of simplicity equilibrium length 
ro was supposed to correspond to the equilibrium value 
of the angle 9 , although it may be not the case. The 
potential W(r) has a form which is typical for non-valent 
interactions, it may be used for description of hydrogen 
bonds and van der Waals interactions. 

A similar model has been used in the papers [22- 
25] for analysis of the longitudinal soliton waves mo- 
tion. It should be mentioned that many polymer macro- 
molecules have 2D zigzag form. For example, in the 
three-dimensional space a polyethylene (PE) molecule 
has a stable plane conformation of trans-zigzag. The 2D 
zigzag form is also usual for hydrogen bond chains ■ ■ -X- 
H- • -X-H- ■ -X-H- • ■ of halides, where X=F, CI, Br, I. 

Let's consider zigzag deformation during its homoge- 
neous stretching. If we denote the extension of the zigzag 
x-step by Ax > 0, then for finding stretched zigzag form 
one has to solve the minimization problem over p 



E{p) = V(p) + U{0) + W(r) -> min : p 



(5) 



with fixed value of x-step h x + Ax [a distance between 
second neighbors is r = 2(h x + Ax), value of the angle 
9 is defined by bond length p\. Solving the problem (5) 
yields the dependence of zigzag energy (per 1 site) from 
the relative x-step extension Ax. 

Let's fix a rigidity of the second neighbors interactions 
K = 1 (interaction energy is eo = k(to — d) 2 /72 = 0.0178), 
while angle deformation energy coefficient e± will be 
changed. 

With ei = the secondary structure of the molecu- 
lar chain is defined only by non-valent interactions be- 
tween second neighbors. Dependence of the angle 9 and 
energy E from the relative uniform extension Ax are 
shown on the Fig. 3. Zigzag stretching is accompanied 
by monotonous angle increasing. When Ax = 0.1836 the 
zigzag straightens out completely in a line (zigzag angle 
is 6 = 180°). 

One can see from the Fig. 3 (b) that the function 
E(Ax) is not convex downward. But if one connects 
two points (Axi, E(Axi)) and (Ax 2 , E(Ax 2 )) (Axi = 
0.0370, Ax 2 = 0.2056) on the curve (Ax, E(Ax)) by the 
straight line which is a tangent of the curve in these two 
points, one gets a convex function. This connecting line 
describes more energetically favorable way of the zigzag 
stretching when one part of the zigzag is in slightly ex- 
tended state with the x-step h x + Axi, and another part 
is in strongly expanded state with the x-step h x + Ax2, 
when the zigzag has straightened out in a line. 

Dependence of the stretching force F = dE/dAx from 
the relative extension Ax is shown on the Fig. 3 (c). 
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Figure 3: Dependence of (a) angle 8, (b) energy E and (c) 
force F = dE/dAx in the zigzag from its relative uniform 
extension Ax. Angle deformation energy coefficient is ei = 
0. The straight line 1 describe a convex hull of the function 
E(Ax). 

One can see from the figure that if non-uniform way of 
stretching is taken into account the force-extension de- 
pendence F(Ax) has a typical form with constant central 
part. When < Ax < Ax± and Ax2 < Ax, that is dur- 
ing uniform stretching of the zigzag, tension increases 
monotonously. When Axi < Ax < Ax%, that is dur- 
ing non-uniform stretching, tension remains constant. A 
typical plateau appears, which has been discovered ear- 
lier in experiments dedicated to stretching of alpha-helix 
[18] and DNA double helix [1, 2]. 

The deformation energy coefficient t\ being increased, 
the typical form of dependencies E(Ax) and F(Ax) is 
preserved up to maximal value ei = 0.015. Here non- 
uniform stretching takes place when 0.0710 < Ax < 
0.1925. Energy benefit of the non-uniform stretching be- 
comes negligible in this case. 

When ei > 0.015 the function E(Ax) becomes convex 
downwards — see the Fig. 4. Here only uniform stretch- 
ing of the zigzag can take place. When t\ = 0.02 zigzag 
angle increases monotonously and under Ax = 0.2343 
reaches the maximal value 8 = 180° (zigzag straightens 
out in a line). The function E(Ax) is convex downwards, 
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Figure 4: The dependence of (a) angle 9, (b) energy E and 
tension F = dE/dAx of the zigzag from a relative extension 
Ax. Angle deformation energy coefficient is e\ — 0.02. 

so the tension F increases monotonously with zigzag ex- 
tension. 

For verifying obtained results let's consider the chang- 
ing of zigzag structure consisting from N = 400 atoms 
while stretching. For this purpose for each value of rela- 
tive extension Ax > one has to solve the minimization 
problem 

N-l N-l N-2 

E v (pJ + E w + E v ^ -> , min „ ( 6 ) 

n=l n=2 n=l W 

with a condition of fixed ends: 

xi = --(N - l)Ax, x N = (N- l)h x + |(JV - l)Ax. 

If {x n ,yn}n=i is a solution of the minimization prob- 
lem (6), then distribution of the longitudinal extension in 
the chain is given by the function Ax n = x n +i — x n — h x , 
and distribution of y-step is given by the function h n = 
\y n +i—yn\- Dependence of distribution Ax n and h n from 
the value of local extension Aa; are shown on the Fig. 5 
and 6. 

One can see from the Fig. 5 that when e\ — a dis- 
tribution of x-step extension Ax n and y-step h n along 
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Figure 5: The dependence of distribution of (a) longitudinal 
extension Ax n and (b) t/-step h n in the stretched zigzag from 
N — 400 atoms with fixed ends from the extension per site 
Ax. Angle deformation energy coefficient is ei = 0. 



the chain corresponds to the behaviour of the function 
E(Ax). When < Ax < 0.037 and Ax > 0.2056 the 
uniform stretching takes place, while when 0.037 < Ax < 
0.2056 the non-uniform one is observed. Here end parts of 
the zigzag are in strongly expanded state with relative ex- 
tension Ax n = 0.2056, while the central part is in slightly 
stretched state with relative extension Ax n = 0.037. On 
the domain boundary gradual transition from one state 
to another is observed. Extension increasing, central 
slightly stretched part monotonously decreases and van- 
ishes with Aa; = 0.2056. Further zizgzag is stretched 
uniformly 

Similar scenario is observed when the limit value t\ = 
0.015 is considered. Here uniform stretching takes place 
when < Ax < 0.0710 and Ax > 0.1925. With 0.0710 < 
Ax < 0.1925 the zigzag has two-phase structure — end 
parts are in strongly expanded state and central part is 
in slightly extended state. 

When e 2 = 0.02, see Fig. 6, the zigzag stretching is 
always uniform. 

So when angle rigidity is zero or rather small stretching 
of a zigzag chain always has the two-phase scenario. As a 
consequence, tension remains constant in some range of 
extensions. Such scenario is due to the typical form of the 




Figure 6: The dependence of distribution of (a) longitudinal 
extension Ax n and (b) j/-step h n in the stretched zigzag from 
N — 400 atoms with fixed ends from the extension per site 
Ax. Angle deformation energy coefficient is e\ — 0.02. 

interaction between second neighbors (4), which have an 
inflection point. So two-phase stretching can be expected 
to occur in molecular chain whose secondary structure is 
formed by weak non-valent interactions. From the other 
hand, if the secondary structure is due to angle and tor- 
sion deformation, then stretching will be uniform, that 
is all sites will be stretched equally. Such a scenario is 
typical for polyethylene (PE) trans-zigzag and polyte- 
trafluoroethylene (PTFE) helix. 

IV. STRETCHING OF THE ALPHA-HELIX 

Let's consider helical molecular chain with atoms hav- 
ing equilibrium coordinates 

R° = (i? cos(n0 o ), Ro sin(n0 o ), nh z ), (7) 

with n = 0,±1,±2,... being atom number, Rq - helical 
radius, 4>o and h z - angle and longitudinal helix period. 
Chain energy is described by the hamiltonian 

H = Y,{l M (Rn^n)+V( Pn )+U(e n )+Z(5 n )+W(r n )}, 

n 

(8) 
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where 3D vector R„(£) stands for the position of the n-th 
site (peptide group) of the helix in the moment t, M is 
site mass. The potential V(p n ) gives an energy of interac- 
tion between neighbor sites n-th and (n + l)-th, depend- 
ing from the distance p n = u n |, where u„ = R n +i — R„ . 
The deformation energy of the angles is described by 
U(6 n ), where 9 n is an angle with a vertex in the n-th 
site, cos(#„) = — (u„, u n _i)//9„_ip„. The third potential 
Z(5 n ) describes deformation energy of the torsion around 
n-th bond. The torsion cosine is given by cos(<5„) = 
i)/l v n|| v n+i|, where the vectorv„ = u n _i x u„. 
The potential W(r n ) stands for the energy of the n-th 
hydrogen bond connecting n-th and (n + 3)-th sites, its 
length is given by r n = |R n+3 - R„|. 

Again for the sake of simplicity we consider dimen- 
sionless model of helix, where equilibrium bond length is 
po = 1, angular step is 0o = 100°, equilibrium angle is 
#o = arccos(-l/3) = 109.47°. Then [26] dimensionless 
helix radius is 

R = po cos(6» /2)/[l - cosOo)] = 0.4919, 

z-stcp is 

Az = pov 7 ! cos(0 o ) + cos(0 o )|/[l - cos(0o)] = 0.6572, 

equilibrium torsion is <5 = arccos(0.2395), equilibrium 
hydrogen bond length is r = 2.0322. 

Bond deformation energy is described by the potential 
(2) with the rigidity K = 10. Angle deformation energy 
is described by the potential (3) with the deformation 
energy coefficient t\ = 1. Torsion deformation energy is 
given by the potential 

W(5) ^[cos^-cos^q)] 2 , (9) 

where £2 > - is deformation energy coefficient. Hy- 
drogen bond energy is given by the potential (4) with a 
parameter d = 0.7 (dimensionless group length) and en- 
ergy eo = k{tq — d) 2 /72 — 0.0246, where rigidity k = 1. 

Scheme of the alpha- helix is shown on the Fig. 7 (a). 
Similar helix models have been used in the papers [26, 27] 
for analysis of ultrasonic soliton motion. 

Let's consider changing of the helix structure while 
stretching. If z-step increases by Az > 0, then to define 
helical structure changing one has to solve a problem of 
finding conditional minimum by two parameters: radius 
R and angular step <p 

E(R, 0, h z + Az) = V{p) + U{6) + W{r) + Z{5) -» min 

(10) 

with the fixed value of z-step h z + Az. Solving the min- 
imization problem (10) yields the dependence of helix 
energy E (per one step) from a value of z-step relative 
extension Az. Firstly the case £2 = is considered. 
Without torsional rigidity helix form is stabilized only 
by hydrogen bonds. Dependencies of the energy E, an- 
gular step <p and radius R of the helix from its relative 




Figure 7: (a) Scheme of the alpha-helix. Helix sites are shown 
on their equilibrium positions, solid lines denote rigid valent 
bonds while dotted lines are for soft hydrogen bonds, (b) A 
zone of the transition between strongly and slightly stretched 
parts of the helix. 



longitudinal extension Az are shown on the Fig. 8. The 
helix being stretched its angular step monotonously in- 
creases and reaches the maximal value cj> = 180° when 
Az = 0.1647 (the helix is unwinded into plane zigzag). 

One can see from the Fig. 8 (a) that the function 
-E(Az) is not convex downward. But if one connects two 
points (Azi, E(Az{)) and (Az 2 , E(Az 2 )) (Azi = 0.0643, 
Az 2 = 0.1694) on the curve (Ax, E(Ax)) by the straight 
line which is a tangent of the curve in these two points, 
it gets a convex function. This connecting line describes 
more energetically favorable way of the zigzag stretching 
where one part of the zigzag is in slightly extended state 
with the z-stcp h z + Azi , and another part is in strongly 
expanded state with the z-step h z + Az2, when the helix 
is unwinded into the plane zigzag. 

When Azi < Az < Az2 the tension of the helix F = 
dE/dAz remains constant. Here the dependence F(Az) 
has typical view with constant central part, as well as for 
the zigzag. 

If torsion deformation energy coefficient e 2 is increased 
this typical view of the dependence E(Az) is preserved 
until the maximum value e 2 — 0.0015 is reached. When 
£2 > 0.0015 the function E(Az) becomes convex down- 
wards - see fig. 9. Here only uniform stretching of the 
helix will take place. When £2 = 0.002 helix stretching 
is accompanied by monotonous increasing of its angular 
step </>, which reaches its maximal value <f> — 180° when 
Az = 0.1727 (the helix has a form of plane zigzag). 

To verify our results about two-phase stretching let's 
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Figure 8: Dependence of (a) the energy E, (b) angular step 
4> and (c) radius R of the helix from its relative longitudinal 
extension Az. Torsion deformation energy coefficient is £2 = 
0. The line 1 shows a convex hull of the function E(Az). 



Figure 9: Dependence of (a) the energy E, (b) angular step 
4> and (c) radius R of the helix from its relative longitudinal 
extension Az. Torsion deformation energy coefficient is £2 = 
0.002. 



consider changing of the structure of stretched helix, con- 
sisting from N = 400 sites, under tension. For each value 
of relative extension Az > one has to solve the mini- 
mization problem 



iV-l 



N-l 



N-2 



N-3 



E v m + E u ( 9 ") + E + E 



71=1 



n=2 



n=2 



71=1 



->■ min : {(R„,, (/>„, z n )}f l=l (11) 
with the condition of fixed ends 

zi =~(N- l)Az, z N = (N- l)h z + i(JV - l)Az. 

The vector (i?„,0„,z„) gives cylindric coordinates of n- 
th site. 

If {(R n , 4>m z n )}n=i i s a solution of minimization prob- 
lem (11), then distribution of longitudinal extension 
along the helix is given by the function Az n = z n+ \ — 
z n ~ h z . Dependence of distribution Az n from relative 
extension Az is shown on the Fig. 10. 



It can be seen from the Fig. 10, that if e 2 = the dis- 
tribution of longitudinal extension Az„ along the chain 
completely corresponds to behaviour of the function 
E(Az). If Az < 0.065 the chain is stretched uniformly, 
while if 0.07 < Az < 0.165 non-uniform stretching is ob- 
served. End parts of the chain are in strongly stretched 
state while the central part is in slightly stretched state. 
On the ends the chain has a form of plane zigzag, in 
the central part it has a form of 3D helix. A bound 
of transition between strongly stretched central part to 
slightly end part is shown on the Fig. 7 (b). Extension 
increasing, central part is monotonously decreased and 
vanishes when Az = 0.170. If extension is larger the 
chain is stretched uniformly. When 62 = 0.0015 non- 
uniform stretching is observed when 0.105 < Az < 0.155. 
If £2 = 0.002 only uniform stretching can be observed. 

Thus independently of torsional rigidity stretching of 
the helix always has two-phase scenario. As a conse- 
quence, tension remains constant in some range of exten- 
sions. Such scenario is due to the typical form of the hy- 
drogen bond potential (4), which have an inflection point. 
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Figure 10: Dependence of distribution of longitudinal exten- 
sion Az n in a stretched helix of N = 400 sites from rela- 
tive extension Az with (a) £2 = 0, (b) £2 = 0.0015 and (c) 
e 2 = 0.002. 

The investigation reveals that extension of polypeptide 
alpha-helix will be non-uniform, what is confirmed by 
all-atomic computer simulations [18]. It should be noted 
that polymer helix which have no hydrogen bonds, such 
as PTFE helix, will be stretched uniformly, tension in- 
creases when extension increases. 



V. STRETCHING OF THE DNA DOUBLE 
HELIX 

For simulating DNA double helix stretching a more 
complex coarse-grained model is used. Each nucleotide 
consists of six grains - three for sugar-phosphate back- 
bone and three for nucleobase. The model is described 
in details in the paper [28]. A scheme of the model is 
depicted on the Fig. 11. In this case a hamiltonian of the 
system has following form: 




Figure 11: Grains involved into valent interactions. Blue lines 
denote valent (harmonic) bonds, red arcs mark valent angles. 
The circles marked as N stand for the atoms N9 on A, G bases 
and N9 on T, C bases (no grain is situated on these atoms, 
their coordinates are calculated additionally from base grains 
positions). 

The first term is a kinetic energy of the system. 

JV 12 ^ 
n— 1 i— 1 

where i is a number of a grain in n-th base pair. 

Terms E v , E a , E t describe energy of correspondingly 
valent bonds, angles and torsions deformation. These 
potentials have rather common form: bond deformation 
energy is calculated by the formula 

U a p(r) = ^K a p(r - R a p) 2 , 

rigidity coefficients K aj 3 and equilibrium values R a p are 
different for different grains. Angle deformation energy 
has a form 

U a (6) = e Q (cos#-cos# Q ) 2 , 

where values of coefficient e a and equilibrium angle differ 
for different angle types. Torsion deformation energy is 
described by the potential 



H — Ek+E v +Eb+E a +E t +E' hb +E st +E e i+E vc iw ■ (12) 



U t = e t (l - cos(0 - O )), 
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where values of coefficient e t and equilibrium torsion an- 
gle 0o differ for different torsions. 

Rotational axis of torsions are shown on the Fig. If. 
The third term Ef, in the hamiltonian (12) describes de- 
formation energy of a nitrogen base. The nitrogen base is 
rather rigid chemical structure what is modeled by rigid 
harmonic potentials which keep four points - grain CI 
and three points on a nucleobase - near one plane. 

The terms E' hb and E st describes interactions between 
nitrogen bases (including stacking and hydrogen bonds) . 
Since a nitrogen base is rather rigid structure, we can 
calculate coordinates of all initial atoms from positions 
of three grains. Thus all-atomic AMBER [29] Coulombic 
and van der Waals potentials are used to mimic hydrogen 
bonds and stacking. 

The last two terms E e i, E v dw describe electrostatic 
and van der Waals interactions between backbone grains. 
Solvent is treated implicitly, using the Generalized Born 
approach [30, 31]. 

Solvent molecules make hydrogen bonds rupture much 
easier (after rupture nitrogen bases form hydrogen bonds 
with water molecules). To take into account this ef- 
fect one has to weaken hydrogen bonds between com- 
plementary bases. For this purposes modified potential 
E' hb = c Q Ehb is used, where E^ is potential of interac- 
tions between complementary bases from all-atomic AM- 
BER potential system [29], and Co = 0.4 is normalizing 
factor. Values of parameters and details of potential cal- 
culation can be found in the paper [28]. 

For definitcness homogeneous sequence polyA-polyT is 
considered. Such double helix can be considered as quasi- 
one-dimensional crystal with the elementary cell being 
one nucleotide pair of the double helix. In the ground 
state each successive nucleotide pair is obtained from its 
predecessor by translation along the z-axis by step Az 
and by rotation around the same axis through helical step 
A(j>. Thus, the energy of the ground state is a function 
of 38 variables: and A(j>, Az, where xxj = 

(xij.i, xij.2, x\.j^) is a radius- vector of j-th grain of the 
first nucleotide pair. 

To find the ground state of the homogeneous double 
helix under tension, it is necessary to solve minimization 
problem over 37 variables 

E = E v + Eb + E a + E t + E' hb + E st + E e i + E v dw 

->min: {xj})^, A<j>, (13) 

under the fixed value of longitudinal step Az (for cal- 
culating the energy one should keep in mind that sum- 
mation is taken over only one base pair and neighboring 
base pairs are obtained from it by rotation and trans- 
lation). As a result, one can obtain the dependence of 
the homogeneous state energy (per one base pair) E on 
the longitudinal step Az - see Fig. 12 (a). One can see 
from the figure that basic minimum is reached with z-step 
Azo = 3.352A, which corresponds to the ground state of 
B conformation of double helix DNA. 

View of ground state under two typical values of z-step 
Az = Azi and Az = Az2 is shown on the Fig. 13. It can 




3 4 5 6 




3 4 5 6 




Figure 12: Dependence of (a) ground state energy E, (b) 
tension F — dE/dAz (s) hydrogen bonds energy —E' hb and 
(d) neighbor basepair stacking interaction energy —E s t from 
longitudinal step Az of polyA DNA ground state. Lines 1,2,3 
give dependencies for uniform stretching while lines 2 and 5 
give dependencies for non-uniform stretching of DNA double 
helix. 



be seen from the figure that in strongly stretched state 
the helix is unwinded more, but all the hydrogen bonds 
remain - only inclination angle of base pair increases. 

One can see from the Fig. 12(a) that when Az > Azq 
the function E(Az) is not convex downward. But if one 
connects two points (Azi, E(Azi)) and (Az 2 , E(Az2)) 
(Azi = 3.75A, Az 2 = 6.15A) on the curve (Ax, E(Ax)) 
by the straight line which is a tangent of the curve in 
these two points, he gets a convex function. This con- 
necting line (Fig. 12, line 2) describes more energetically 
favorable way of the helix stretching where one part of 
the zigzag is in slightly extended state with the z-step 
Azi , and another part is in strongly expanded state with 
the z-step Az 2 . The fig. 12 shows that such scenario 
consumes less energy than uniform stretching of all sites. 

Dependence of hydrogen bond energy from extension is 
shown on the Fig. 12 (c). Under tension hydrogen bonds 
weaken. When z-step Az = Azi hydrogen bond weakens 
approximately twice. One can see from the Fig. 12 (d) 
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Figure 13: View of ground state of (a) slightly stretched (lon- 
gitudinal step is Az = Azi) and (b) strongly stretched (z-step 
is Az = Az 2 ) polyA DNA double helix. 




Figure 14: View of ground state of (a) slightly stretched (lon- 
gitudinal step is Az = Azi) and (b) strongly stretched (z-step 
is Az — AZ2) polyA DNA double helix under the temperature 
T = 300 K. 



that weakening of stacking interactions also occurs while 
stretching. Thus we can conclude that both hydrogen 
bonds and stacking interactions are responsible for non- 
convexity of the function E(Az). 

Under the room temperature T = 300K weakened hy- 
drogen bonds rupture due to thermal fluctuation - see 
Fig. 14. Therefore in the case of non-uniform stretching 
DNA macromolcculc consists of two parts: one part is 



slightly stretched double helix where hydrogen bonds re- 
main [see Fig. 14 (a)], another one is strongly stretched 
double helix where some hydrogen bonds are ruptured 
[see Fig. 14 (b)] . Exactly this behaviour is observed in all- 
atomic DNA stretching simulations[15, 32]. In the papers 
[3, 4, 9, 10] existing of the plateau in force-extension dia- 
gram is attributed to the rupture of hydrogen bonds. But 
our analysis reveals that it isn't the case. The plateau 
is observed also in the absence of thermal fluctuations, 
when hydrogen bonds are only weakened, but not yet 
ruptured. 

Force-extension dependence F = dE(Az) /dAz is 
shown on the Fig. 12 (b). The force is seen to increase 
monotonously in the range of z-step from Azq to Azi 
(curve 3). When the value of z-step Azi is reached uni- 
form stretching became non-uniform one. In the range 
Azi < Az < Az2, when non- uniform stretching takes 
place, tension remains constant (curve 5). When the 
value Az2 is reached stretching again becomes uniform, 
and tension increases monotonously with Az increasing 
(curve 4). 

The plateau in force-extension diagram has been ob- 
served previously in single molecule stretching experi- 
ments [1, 2], the plateau was reported to occur in the 
range 1.1 < Az/Az < 1.7. So, having Az = 3.34, 
experiments estimate plateau range as 3.67A< Az < 
5.68A, which agree well with current results, obtained 
within the framework of CG representation 3.75A< Az < 
6.15A [see Fig. 12 (b)]. Value of there transition force 
0.2 eV/A= 320 pN is somewhat higher than experimen- 
tal estimates, this matter is discussed below. It should 
be noted that this value should be sequence-depended. 

Molecular dynamics simulation have demonstrated 
that the scenario of non-uniform stretching is preserved 
under the room temperature 300K, only range of the 
plateau increases slightly, but the value of force coincide 
with the value under T=0 K. 

In the reported experiments the value of tension in 
the plateau depends significantly from conditions, such 
as temperature, salt, sequence and especially from a way 
of fixing ends. The typical value of force is usually re- 
ported to be about 65 -j- 70 pN [1, 2, 12]; however one 
should keep in mind that this value was obtained while 
stretching by the same type (3' or 5') ends, while the 
other two ends were free. On the contrary, the simu- 
lations reported here corresponds to uniform stretching 
by all four ends. Experimentally, larger values of the 
plateau tension were reported in this case: 105 120pN 
[12, 33], although that is is still less than half of the 320 
pN predicted by our model. However, since many DNA 
properties are strongly sequence-dependent, only semi- 
quantitative comparison of our homogeneous model with 
the above experiment (that used a specific, non-uniform 
sequence) is possible. It is known from other experiments 
that, for example, pure poly(dG-dC) and poly(dA-dT) 
DNA sequences, pulled by two of the same type ends, 
yield tension values at the plateau that differ by a factor 
of two [33]. 
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It should be mentioned that difference between results 
of all-atomic simulations and experimental values is of 
the same order [1, 13]. Therefore numerical agreement 
with experimental values can be considered as reason- 
able. 



VI. CONCLUSIONS 

The investigation demonstrates that existing of non- 
valent interactions can result in appearing an inflection 
point in the dependence of homogeneous chain energy 
from extension. As a consequence a non-uniform way 
of stretching becomes energetically favorable, when one 
part of a chain is stretched slightly and another part is 
stretched strongly. Extension of the whole chain is due 
to increasing of fraction of strongly stretched sites, so 
the tension remains constant. Such way of non-uniform 
stretching explains the plateau in a force-extension dia- 
gram. 

The non-uniform stretching is observed in stretching 
simulations of both polypeptide alpha-helix and DNA 



double helix. Obtained results enable to give a general 
explanation of the force-extension plateau observed both 
for DNA double helix [1, 2, 13] and polypeptide alpha- 
helix [6, 7, 18] in single molecule stretching experiments, 
as well as in corresponding computer simulations. In the 
case of alpha-helix gradual unwinding of the helix have 
been supposed to occur while stretching [6]. However 
suggested models of this event were based upon the as- 
sumption that strongly stretched state is stable, that is 
corresponds to the free energy minimum. Our investiga- 
tions reveal that it may be not the case. 

In the case of DNA molecule the force-extension 
plateau is explained either by cooperative transition into 
the S-form conformation [1, 13, 14] or by force- induced 
melting [3, 4, 9, 10]. Our results allow to conclude that 
stretching is accompanied by gradual transition from 
slightly stretched to strongly stretched state, but the 
strongly stretched state isn't stable. However some hy- 
drogen bonds ruptures may occur due to thermal fluc- 
tuations. Obtained estimates of force-extension plateau 
are in reasonable agreement with experiments. 
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